arXiv:1507.06238vl [cond-mat.stat-mech] 22Jul2015 


The Maximum of a Fractional Brownian Motion: Analytic Results from Perturbation 

Theory 

Mathieu Delorme and Kay Jorg Wiese 

CNRS-Laboratoire de Physique Theorique de I’Ecole Normale Superieure, 24 rue Lhomond, 75005 Paris, France 


Fractional Brownian motion is a non-Markovian Gaussian process Xt, indexed by the Hurst 
exponent H. It generalises standard Brownian motion (corresponding to H = 1/2). We study the 
probability distribution of the maximum m of the process and the time t m ax at which the maximum is 
reached. They are encoded in a path integral, which we evaluate perturbatively around a Brownian, 
setting H = 1/2+e. This allows us to derive analytic results beyond the scaling exponents. Extensive 
numerical simulations for different values of H test these analytical predictions and show excellent 
agreement, even for large e. 


Random processes are ubiquitous in nature. While av¬ 
eraged quantities have been studied extensively and are 
well characterized, it is often more important to under¬ 
stand the extremal behavior of these processes [I] , asso¬ 
ciated with failure in fracture or earthquakes, a crash in 
the stock market, the breakage of dams, etc. 

Though many processes can successfully be modeled 
by Markov chains and are well analyzed by tools of sta¬ 
tistical mechanics, there are also interesting and realis¬ 
tic systems which do not evolve with independent incre¬ 
ments, and thus are non-Markovian, i.e. history depen¬ 
dent. Dropping the Markov property, but demanding 
that a continuous process be scale-invariant and Gaus¬ 
sian with stationary increments defines an enlarged class 
of random processes, known as fractional Brownian mo¬ 
tion (fBm). Such processes appear in a broad range of 
contexts: Anomalous diffusion [2], polymer translocation 
through a pore the dynamics of a tagged monomer 

0 0 , finance (fractional Black-Scholes and fractional 
stochastic volatility models 0 ), hydrology [9], and many 
more. 

FBm is a generalization of standard Brownian mo¬ 
tion to other fractal dimensions, introduced in its final 
form by Mandelbrot and Van Ness El- It is a Gaussian 
process (X t )te r> starting at zero, Vo = 0, with mean 
(Xt) = 0 and covariance function (variance) 

(X t X s )=s 2H +t 2H -\t-s\ 2H . (1) 

The parameter H £ (0,1) is the Hurst exponent; the 
process typically grows with time as t H . Standard Brow¬ 
nian motion corresponds to H = 1/2; there the covari¬ 
ance function reduces to (X f X s ) = 2min(s,f). Unless 
H = 1/2, the process is non-Markovian, i.e. its incre¬ 
ments are not independent: For H > 1/2 they are corre¬ 
lated, whereas for H <1/2 they are anti-correlated, 

(d t X t d s X s ) = 2H(2H - l)|s - . (2) 
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FIG. 1: Two realisations of fBm paths for different values of 
H, generated using the same random numbers for the Fourier 
modes in the Davis and Harte procedure m- The observables 
m and t max are given. 


X 0 = 0 and total time T > 0. Figure [l] shows an illus¬ 
tration for different values of H , using the same random 
numbers for the Fourier modes. We will denote Pj)(m) 
and Pjj(t) their respective probability distributions. Pre¬ 
vious studies can be found in mm- 

These observables are closely linked to other quanti¬ 
ties of interest, such as the first-return time, the survival 
probability, the persistence exponent, and the statistics 
of records. Though studied since long time, most results 
for non-Markovian processes are quite recent mm- 
Following the ideas of we encode our observ¬ 

ables Pjj(m) and Pjj(t) in a path-integral, 

Z + (m 1 ,ti\x 0 -,m 2 ,t 2 ) = 

px tl +t 2 =rn 2 ( 3 ) 

/ T>[X]0[X]6(X tl -x o)e- s A]. 

J Xq—TTIi 


In this letter we study the maximum of a fractional Brow¬ 
nian motion m = max t6 [ 0iT ] X t and the time f max when 
this maximum is reached [2T] with the initial condition 


This sums over all paths X t , weighted by their prob¬ 
ability e~ s ^ x \ starting at X 0 = m\ > 0 (shifted for 
convenience), passing through Xq (close to 0) at time 
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FIG. 2: Graphical representation of a contribution to the 
path-integral Z + (mi, t\\ xo', m 2 , £ 2 ) given in Eq. ©• The red 
curve represents the non-local interaction in the action (sec¬ 
ond line of Eq. 0 ) while blue lines are bare propagators. 
There are two other contributions when the time ordering is 
ti < T 2 < fi or fi < n < T 2 , already computed in m- 


p T H {t) 
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FIG. 3: Distribution of f max for T = 1 and H = 0.25 (red) 
or H = 0.75 (blue) given in Eq. dTp (plain lines) compared to 
the scaling ansatz, i.e. T = cst (dashed lines) and numerical 
simulations (dots). For H < 0.5 realisations with t max ~T/2 
are less probable (by about 10%) than expected from scaling. 
For H > 0.5 the correction has the opposite sign. 


t\. and ending in X tl+t2 = m 2 > 0, while staying 
positive for all t £ [0,fi + f 2 ]■ The latter is enforced by 
the product of Heaviside functions 9[X} := JlltJo* 2 0(X S ). 


As X t is a Gaussian process, the action S can (at least 
formally) be constructed from the covariance function of 
X t . However this is not enough to evaluate the path 
integral © in all generality. Following the formalism of 
nn, we use standard Brownian motion as a starting point 
for a perturbative expansion, setting H = f + e with e 
a small parameter; then the action at first order in e is 
(we refer to the appendix of [19] for the derivation) 


S[X\ 
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© is linked to this distribution via 
Pjj(t) = lim [ Z + (mi,t-x 0 -m 2 ,T-t). (6) 

x 0 ->0 Zi Jmi,m 2 >0 


The normalization Z depends on Xq and T. Our result 
for the distribution of f max takes a nice forme if we ex¬ 
ponentiate the order-e correction obtained from Eq. 

P « (t)= [KT-Dp e * p ( e ' F (^)) +0(e2) • (7 > 

This is plotted on Fig. [3] We see the expected change in 
the scaling form of the Arcsine law, \Jt(T — t) —> [t(T — 
t)] H and a non-trivial change in the shape given by the 
function 


F{u) =\fu [7T — 2 arctan(v / rt)] 


7r — 2 arctan 


+ cst . 


( 8 ) 


The time r is a regularization cutoff for coinciding times 
(one can also introduce discrete times spaced by r iisi)- 
The first line is the action for standard Brownian motion, 
with a rescaled diffusion constant [25] D e T = 1 + 2e(l + 
ln(r)) + 0(e 2 ) ~ (er) 2e . The second line is a correction, 
non-local in time since fBm is non-Markovian. 

Computing the e expansion of ^ using Q is rather 
technical. A graphical representation of the key term 
is given in Fig. [2J The result for Z + (m±, fy; Xo; to 2 , f 2 ) 
covers a page, presented in [20]. We use this result here 
to deduce its most interesting implications, starting with 
the probability distribution of t = f max . For Browanian 
motion (H = 1/2), this distribution is well known as the 
Arcsine law [26] . 

b rw = 75(bir for<eM ' (5) 

Until now, only scaling properties were known for this 
distribution in the general case m- The path integral 


The time reversal symmetry t —> T — t (corresponding to 
u —)■ ia^ 1 ) is explicit; the constant ensures normalization. 

We tested the prediction ©-© with numerical sim¬ 
ulations of a discretized fractional Brownian motion for 
different values of H. To this aim, we used the Davis 
and Harte procedure as described in [JO] (and references 
therein). To compare numerical results with the theory, 
we extract an estimation of the function J- as 

^um := \ !n (pLI (i) X [t(T - i)]*) . (9) 

Here ///;// (t) is the numerical estimation of the distribu¬ 
tion of t max for the discretized fBm at given H (obtained 
with uniform binning). Apart from discretization effects, 
we should see significant statistical errors as e —► 0, and 
systematic order-£ 2 corrections for larger e. As can be 
seen on Figs. [3] and [4] our numerical and analytical re¬ 
sults are in remarkable agreement for all values of H stud¬ 
ied, both for e positive and negative. As an example, for 
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FIG. 4: Left: Numerical estimation of T for different values of H on a discrete system of size N = 2 12 , using 10 s realizations. 
Plain curves represent the theoretical prediction |8|, vertically translated for better visualization. Error bars are 2<r estimates. 
Note that for H = 0.6, H = 0.66 and H = 0.8 the expansion parameter e is positive, while for H = 0.4, H = 0.33 and H = 0.2 
it is negative. Right: Deviation for large |e| between the theoretical prediction (|8]l and the numerical estimations (|9|, rescaled 
by e. These curves collapse for different values of H , allowing for an estimate of the 0(e 2 ) correction to 


H = 0.75, the correction to the pure scaling distribution 
has a relative magnitude of 10% (see Fig. [3]), which is 
measured in our simulation with a relative precision of 
0.5%. This precision even allows to numerically extract 
the subleading 0(e 2 ) correction, see Fig. fright. 

We now present results for the distribution of the max¬ 
imum For Brownian motion 

m 2 
P 4T 

P I M = —j= , m > 0 . (10) 

On the other hand, not much is known for generic values 
of H. This distribution is of interest, as it is linked to 
the survival probability S(T,x), and the persistence ex¬ 
ponent 9. The latter is defined for any random process 
X t with Xq = x as 

S(T, x) = proba (X t > 0 for all t £ [0, T ]) 

T -o. _ (11) 

T —>-oo 


For a large class of processes the exponent 9 is indepen¬ 
dent of x. For fractional Brownian motion with Hurst 
exponent H it was shown that 9 X = 9 = 1 — H mm- 
To understand the link of S(T , x) with the maximum dis¬ 
tribution for fBm, we use self affinity of the process X t 
to write P^(m) as 


p h M = 


y/2T 


H 


fn 


m 


V2T 


H 


( 12 ) 


Here / is a scaling function depending on H. Eq. (10) 
can be reformulated as fi(y) = Y^e -y2 / 2 . The survival 
probability is related to the maximum distribution by 


r x r J. c t h 

S(T,x)= / P T (m)din= / fn{u)du . 
Jo Jo 


(13) 


This states that a realisation of a fBm starting at x and 
remaining positive is the same as a realisation starting at 
0 with a minimum larger than —x, due to translation in¬ 
variance of the fBm. Finally, the symmetry x —> —x gives 
the correspondence between minimum and maximum. 

These considerations allow us to predict the scaling 
behavior of P^ ( m ) at small m from the large-T behaviour 
of S(T, x) [15] . 

Pft(m) ~ mTi- 1 =mTt~ 2 . (14) 


Using our path integral, we can go further. The maxi¬ 
mum distribution can be extracted from Eq. ([3]), 


P]j(m) = lim — ( dt f Z + (m, t; Xo; m2, T — t) . 
x-o-m Z J Q Jm 2 >0 

(15) 

Its e-expension leads to the scaling form of Eq. (12), with 


fH(y) 


e -4- e e[0(»)+cst] 


+ 0(e 2 ). 


(16) 


The constant term ensures normalization. The function 
Q involves the hypergeometric function 2 F 2 : 

5( ») 4,4,44)-v 

+ tt(1 - y 2 ) erfi. 

+ {y 2 - 2 ) [ 7E + ln(2j/ 2 )] ■ 

This function has a different asymptotics for small and 
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FIG. 5: Left: The combination (21 1 for H = 0.6. The plain line is the analytical prediction exp(e[5(m/v / 2) + 41nm] + cst) 
of the distribution of the maximum without its small-scale power law and large-scale Gaussian behavior. The symbols are 
numerical estimations for T = 1 of the same quantity m 2 ^ 1 ' H exp(m 2 /(m) for different sample sizes. At small scale 
discretization errors appear. At large scales the statistics is poor due to the Gaussian prefactor. For the four decades in between 
theory and numerics are in very good agreement. Right: ibid for H = 0.75. The large scale-behavior on both plots is consistent 
with m 2e . 


large y, 


G(y) ~ - 2 ln(y) 

y—> oo 

~ -41n(y) . 
v-y o 


(18) 


The second line implies that Pjj(m) ~ m 4e which 


m—> 0 


is consistent (at order e) with the scaling result (141, 


jj— 2 = —4 £+0(s 2 ). Formulas (16)-(17) also predict the 
distribution at large m. The leading behavior of P^(m) 
is Gaussian, which is well known, and can be derived from 
the Borrel inequality [23]. Our result for the subleading 
term can be written as 


lim 

y—too 


In (/i?(j/)exp(^)) 


My) 


= -2 e + 0(e 2 ). (19) 


In order to test these predictions against numerical sim¬ 


ulations, we can rewrite the form (16) s.t. the small-m 


behavior matches the exact scaling result (14) 


f H (y) = ^ e -^yTt- 2 e E ^ +Alnv+cst] +0{e 2 ) . 


( 20 ) 


To extract the non-trivial contribution from numerical 
simulations, we study for T = 1 (see Fig. [5]) 

m 2 -^e^P n ^(m) =e £ [ e (^) +41nm+cs t] + 0(e 2 ) . 

_ ( 21 ) 

The sample size N (be. lattice spacing dt = N 1 ) of 
the discretized fBm used for this numerical test is im¬ 
portant, as Pium ( m ) recovers Brownian behavior for m 
smaller than a cutoff of order N~ H . Far small H the nec¬ 
essary system size is very large, so we focus on H > 0.5. 


Figure [5] presents results for H = 0.6 and H = 0.75, 
without any fitting parameter. The constant term in the 
scaling form, relevant for normalization, is evaluated nu¬ 
merically. As predicted, convergence to the small-scale 
behavior is quite slow. This would lead to a wrong nu¬ 
merical estimation of the persistence exponent or other 
related quantities if the crossover to the large-scale be¬ 
havior is not properly taken into account. At large scales, 
the numerical data on Fig.[5]grow as m 2e , consistent with 
the prediction (19). 

To conclude, we have given analytical results for the 
maximum of a fractional Brownian motion, and the time 
when this maximum is reached. To our knowledge these 
are the first analytical results for generic values of H in 
the range 0 < H < 1, beyond scaling relations. Compar¬ 
ison to numerical simulations shows excellent agreement, 
even far from the expansion point H = 

Our calculations also gave the joint probability of the 
maximum, the time when the maximum is reached, and 
the final point [20:. This allows us to address other 
observables of interests, such as fractional Brownian 
bridges. 
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